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' Abstract 

A fast algorithm for Antoine and Vandergheynst's (1998) directional Continuous Spherical 
Wavelet Transform (CSWT) is presented. Computational requirements are reduced by a 
> ■ factor of O(-yjNpix), when A p ; x is the number of pixels on the sphere. The spherical Mexican 

Hat wavelet Gaussianity analysis of the WMAP 1-year data performed by Vielva et al. (2003) 
is reproduced and confirmed using the fast CSWT. The proposed extension to directional 
analysis is inherently afforded by the fast CSWT algorithm. 




1 Introduction 



A range of primordial processes imprint signatures on the temperature fluctuations of the Cosmic 
Microwave Backround (CMB). For instance, various cosmic defect and non-standard inflationary 
models predict non-Gaussian anisotropics. By studying the Gaussianity of the CMB anisotropics 
evidence may be provided for competing scenarios of the early Universe. In addition, a num- 
ber of astrophysical processes introduce secondary sources of non- Gaussianity. Measurement 
systematics or contamination may also be highlighted by Gaussianity analysis. 

Wavelets are a powerful tool for probing the Gaussianity of CMB anisotropies. Previous 
wavelet analysis of the CMB, however, has been restricted to simple spherical Haar and isotropic 
Mexican Hat wavelets. A directional analysis on the full sky has previously been prohibited by 
the computational infeasibility of any implementation. We rectify this problem by providing 
a fast algorithm for Antoine and Vandergheynst'JP Continuous Spherical Wavelet Transform 
(CSWT), based on the fast spherical convolution proposed by Wandelt and GorskP. 

The remainder of this paper is organised as follows. The CSWT is presented in Section [21 
and the fast implementation in Section |31 In Section 0] the fast CSWT is applied to reproduce 
the non-Gaussianity CMB analysis performed by Vielva et aJ3. Concluding remarks are made 
in Section [5J 



2 A directional Continuous Spherical Wavelet Transform 



Antoine and VandergheynspJ extend Euclidean wavelet analysis to spherical geometry by con- 
structing a wavelet basis on the sphere. The natural extension of Euclidean motions on the 
sphere are rotations, defined by (JZ p f)(uj) = f(p~ l uj), p £ 50(3), where we parameterise p by 
the Euler angles (a, {3, 7). Dilations on the sphere, denoted (V a f)(u) = f a (io), are constructed 
by first lifting the sphere 5 2 to the plane by a norm preserving stereographic projection from the 
South pole, performing the usual Euclidean dilation in the plane, before re-projecting back onto 
5 2 . Mother spherical wavelets are constructed simply by projecting Euclidean planar wavelets 
onto the sphere by a norm preserving inverse stereographic projection. A wavelet basis on 5 2 
may be constructed from rotations and dilations of an admissible mother spherical wavelet. The 
corresponding wavelet family {i^ a ,p = T^pT^aip, P S 50(3), a G K+} provides an over-complete 
set of functions in L 2 (5 2 ). The CSWT is given by the projection onto each wavelet basis function 

S(a,a,/3,j)= (7£q, fl 7 V> a )*(w) s(u) dp(ui) , (1) 

where the * denotes complex conjugation and dp{tu) = sin(0) d9 d<p is the usual rotation invariant 
measure on the sphere. 

3 Fast algorithm 

A direct implementation of the CSWT is simply not computationally feasible for a data set 
of any practical size; a fast algorithm is essential. At a particular scale the CSWT is essen- 
tially a spherical convolution, hence it is possible to apply Wandelt and Gorski'0 fast spherical 
convolution algorithm to rapidly evaluate the transform. 

3.1 Fast implementation 

There does not exist any finite point set on the sphere that is invariant under rotations, hence 
it is more natural, and in fact more accurate for numerical purposes, to recast the CSWT in 
spherical harmonic space. The Wigner rotation matrices (defined by Brink and SatchlePJ for 
example) introduced to characterise the rotation of a spherical harmonic may be decomposed 
as D l mm ,(a, f3,j) = e~ ima d l mm ,(f3) e~ tmi . This decomposition is exploited by factoring the 
rotation into two separate rotations, both of which contain a constant ±7r/2 polar rotation: 
tt/2, -tt/2, (3 7^-0, 71-/2, 7+7r/2- By uniformly sampling and applying some algebra the 
CSWT may be recast as 

j=0 j'=o j"=0 

where the summation is simply the unnormalised 3D inverse discrete Fourier transform of 

^max 

t m+l max ,m'+l max ,m"+m max = e l{m ~ m)7r/2 ^ d^, m (7r/2) 4n'm" C 71 "/ 2 ) ^n" * > (3) 

Z=max([m|,|m'|,|m"[) 

where^ m denote spherical harmonic coefficients, l ma x and rn max define the general and azimuthal 
band limits of the wavelet respectively and the shifted indices show the conversion between the 
harmonic and Fourier conventions. The CSWT may be performed very rapidly in spherical 
harmonic space by using fast Fourier techniques to rapidly and simultaneously evaluate (J2J), 
once © is precomputed. a 

"Memory and computational requirements may be reduced by a further factor of two for real signals by 
exploiting the conjugate symmetry relationship f_ m _ m / _ m » = t* m m , m n. 



Table 1: CSWT algorithm comparisons 



Algorithm 


Complexity 


Direct 

Semi-fast 

Fast 


C(L"iV 7 ) 
0(L 3 log 2 (L)iV 7 ) 
G>(L 3 iV 7 ) 



(a) Complexity: Let the number of samples 
of each parameter be of order L, except iV 7 
which typically may be much lower. 



N pix 


Execution time 

(mimsec) 

Direct Semi-fast Fast 


768 
3,072 
12,288 
786,432 


00:01.19 00 
00:18.60 00 
05:01.48 04 


01.12 00:00.01 
17.38 00:00.04 
43.06 00:00.21 
01:54.15 



(b) Typical execution times: Tests performed on a 
dual 900MHz processor system with 4GB of memory. 



3.2 Comparison with other algorithms 

Direct and semi-fast (where only one rotation is performed in Fourier space) implementations of 
the CSWT are also possible. A comparison of the theoretical complexity and typical execution 
times of each algorithm is presented in Tabled The fast CSWT algorithm provides a saving of 
0(^N p i x ) for iVpi x pixels on the sphere. 

4 CMB non-Gaussianity analysis 

We reproduce the Gaussianity analysis of Vielva et aJ3, preprocessing the WMAP data in the 
same manner. The resolution of the co- added map defined by Komatsu et cJ3 is down-sampled 
by a factor or 4, before the KpO exclusion mask is applied to remove emissions due to the 
Galactic plane and known point sources. 

4-1 Spherical wavelet analysis 

The CSWT is a linear operation, hence the wavelet coefficients of a Gaussian map will also 
obey a Gaussian distribution. To test for deviations from Gaussianity, skewness and kurtosis 
statistics are calculated for each wavelet coefficient map at each scale. Monte Carlo simulations 
are performed to construct confidence bounds on the test statistics. 

The application of the KpO mask distorts coefficients corresponding to wavelets that overlap 
with the mask exclusion region. These wavelet coefficients must be removed from any subse- 
quent analysis. Our construction of an extended coefficient exclusion mask differs to that of 
Vielva et alP and inherently accounts for the dominant distortion (either point-source or Galac- 
tic plane) at each scale. The only non-zero coefficients in a CSWT of the original mask are those 
that are distorted (due to the zero-mean property of spherical wavelets). These may be easily 
detected and the coefficient exclusion mask extended accordingly. 

4-2 Results 

We reproduce the results of Vielva et aJ3 for the spherical Mexican Hat wavelet analysis of the 
co-added WMAP data. The Mexican Hat wavelet scales {ai}]^ = {14,25,50,75,100,150,200, 
250, 300, 400, 500} ji x acrmin are considered, corresponding to an effective size of the sky of 
d = 4tan~ 1 (ai/v / 2) 2y2aj (defined as the angular separation between opposite zero-crossings). 
Figure ^ shows the skewness and kurtosis of the coefficients at each scale. The wavelet analysis 
inherently allows one to localise signal components on the sky, as illustrated in Figure EJ We 
make similar observations to Vielva et aJ3, although the different coefficient exclusion masks 
produce slight discrepancies. These discrepancies do not alter the general findings of the analysis. 



(a) Skewness 



(b) Kurtosis 



Figure 1: Spherical Mexican Hat wavelet coefficient statistics: Confidence regions derived from Monte Carlo 
simulations are shown for 68% (red), 95% (orange) and 99% (yellow) levels, as is the mean (solid white line). 




(a) Original coefficient map (b) Thresholded coefficient map 

Figure 2: Spherical Mexican Hat wavelet coefficients at scale ag = 250': Those coefficients below 
3a(as) are thresholded to zero so that likely deviations from Gaussianity may be localised on the CMB sky. 

5 Conclusions and future work 

A fast algorithm is presented and evaluated for performing a directional CSWT on the sphere. 
The fast implementation reduces the complexity of the CSWT by 0{^J N p i^) : where iVpi x is the 
number of pixels on the sphere. Furthermore, the numerical accuracy of the CSWT is improved 
by elegantly representing rotations in harmonic space. 

The Gaussianity analysis of the WMAP 1-year data performed by Vielva et alP has been 
reproduced and confirmed using the fast CSWT. We consider the extension to a full directional 
analysis in an upcoming publication by McEwen et preliminary findings indicate deviations 
from Gaussianity outside of the 99% confidence level. 
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